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A particular periodic orbit in the Earth-Sun circular restricted three body problem is 
shown to have the characteristics needed for a ballistic lunar capture transfer. An injection 
from a circular parking orbit into the periodic orbit serves as an initial guess for a targeting 
algorithm. By targeting appropriate parameters incrementally in increasingly complicated 
force models and using precise derivatives calculated from the state transition matrix, a 
reliable algorithm is produced. Ballistic lunar capture trajectories in restricted four body 
systems are shown to be able to be produced in a systematic way. 
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Nomenclature 

weak stability boundary 
low Earth orbit 

circular restricted three body problem 

restricted three body problem 

restricted four body problem 

gravitational parameter of the Sun 

gravitational parameter of the Earth 

gravitational parameter of the Moon 

mass ratio in the circular restricted three body problem 

Jacobi constant in the circular restricted three body problem 

orientation angle of a circular parking orbit 

orientation angle of a circular parking orbit 

orientation angle of a circular parking orbit 

magnitude of the transfer trajectory insertion burn 

position vector of a spacecraft in a coordinate system centered on the Earth-Moon barycenter 
velocity vector of a spacecraft in a coordinate system centered on the Earth-Moon barycenter 
acceleration vector of a spacecraft in a coordinate system centered on the Earth-Moon barycenter 
position vector of the Moon in a coordinate system centered on the Earth-Moon barycenter 
position vector of the Sun in a coordinate system centered on the Earth-Moon barycenter 
position vector of L 2 

position vector of a spacecraft relative to the Sun 
position vector of a spacecraft relative to the Earth 
position vector of a spacecraft relative to the Moon 
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r, = distance of L-> from the Moon 

‘-‘2! m 

x p = vector of free parameters in the targeting and optimization algorithms 

c = constraint vector in the targeting and optimization algorithms 

J = performance index in the optimization algorithm 

KE m = Keplerian energy with respect to the Moon 

X = state vector 

<D(t f ,t 0 ) = state transition matrix 

F(t) = state propagation matrix 

f m = true anomaly in a Moon-centered coordinate frame 

Introduction 

A ballistic lunar capture trajectory is characterized by a spacecraft transitioning from a hyperbolic lunar orbit into 
an elliptical lunar orbit without the need of a maneuver at the time of transition. The trajectories were initially 
developed using the concept of the weak stability boundary [1, 2] (WSB), a chaotic region of phase space that exists 
around masses in multi-body problems. Trajectories of this type were first developed by Belbmno in 1987 [1]. 
Belbmno and Miller [2] were able to use the WSB in combination with the Sun’s gravitational influence to 
successfully place a satellite in lunar orbit in 1991. Presented here is a new method of targeting impulsive ballistic 
lunar capture transfers from low Earth orbit using a particular family of periodic orbits in the Sun-Earth circular 
restricted three body problem (CRTBP). 

Several methods of designing ballistic lunar capture transfers have been documented [3-11]. Belbmno and 
Carrico [3] have presented a two dimensional, user-in-the-loop forward targeting algorithm for transfers that 
approach the Moon along an orbit that has an apogee that is much greater than the distance from the Earth to the 
Moon. Yamakawa, et al. [4] used a shooting method to produce similar transfers in the planar restricted four body 
problem (RFBP). Additional algorithms exist that utilize the invariant manifolds associated with libration point 
orbits in the Sun-Earth and Earth-Moon circular restricted three body problems. For example, Koon, et al. [5] 
targeted ballistic lunar capture trajectories by finding intersections of these invariant manifolds. Similarly, Yamato 
and Spencer [6] approximated the invariant manifolds in a perturbed CRTBP, yielding transit orbits. Parker and Lo 

[7] have categorized families of ballistic lunar capture trajectories found from invariant manifolds with the intention 
of allowing a mission planner to choose appropriate trajectories for specific missions. Alternatively, Biesbroek et al. 

[8] have used genetic algorithms to successfully find WSB trajectories. Finally, Yagasaki [9] has created a non- 
linear boundary value problem that obtains a solution by beginning with an elliptic arc in the two body problem 
(mass parameters of the Sun and Moon equal to zero), and iterating the solution with increasing mass parameters of 
the Sun and Moon until the real-world solution is obtained. 



Other low energy targeting methods exist that do not rely on solar gravitational effects. The trajectories that 
these transfers target approach the Moon from the direction of the interior Earth-Moon Lagrange point. Bolt and 
Meiss [10] developed a targeting scheme in the Earth-Moon planar circular restricted three body problem that relies 
on recurrence of chaotic trajectories. Macau and Grebogi [1 1] used a similar method to target transfers to the moon 
through chaotic spaces in the restricted three body problem through elimination of recurrent orbits, however, 
improved transfer times over Bolt and Meiss were achieved at the expense of a second maneuver and higher fuel 
costs. Both of these methods rely on very large parking orbits around the Earth to achieve a low energy impulsive 
transfer without the need of a solar perturbation. Mengali and Quarta [12] also disregard the solar gravitational 
influence when they compare their planar three body bi-impulsive method to WSB transfers. 

The targeting of ballistic lunar capture trajectories poses problems because of the chaotic nature of the trajectory 
and the complexity of the four body dynamics that are required to produce them. The WSB is by definition a chaotic 
region in phase space where small changes in the state of a spacecraft will lead to large variations in its trajectory. In 
addition, ballistic lunar captures that rely on a single main impulsive maneuver to transfer from the Earth to the 
Moon require the gravitational dynamics of the Sun, Earth, and Moon. The four body dynamical system makes 
systematic convergence to ballistic lunar capture trajectories problematic. 

Attempts have been made to use simpler dynamics to understand ballistic lunar capture trajectories for the 
purpose of creating initial guesses accurate enough to converge to solutions in a numerical targeting scheme. 
Ivashkin, [13] for example, has compared them to bi-elliptic transfers in a central body gravity field. Here, the solar 
perturbation of the Earth-centered two body problem provides the intermediate AV that increases the perigee 
distance of the orbit. From Lidov, [14] an estimation of the effect of this perturbation on the perigee distance can be 
made. Although the comparison of a bi-elliptic transfer with a ballistic lunar capture trajectory has been shown by 
these authors to be geometrically appropriate, the dynamics of the ballistic lunar capture trajectory can be better 
demonstrated in a three body model. 

The solar effect on perigee can be seen in a family of periodic orbits in the Sun-Earth circular restricted three 
body problem. These periodic orbits explicitly exhibit the influence of the Sun on the low energy transfers. A 
targeting algorithm is described that uses a particular member of a family, documented by Markellos [15] as family 
f\ 6, as a generating trajectory for an Earth-Moon transfer in the restricted four body problem. The algorithm 
overcomes the stated difficulties in numerically seeking ballistic lunar capture trajectories by using an incremental 



approach from an appropriate initial guess and precise derivatives based on the state transition matrix. The given 
parameters for the algorithm are the initial transfer date and the properties of a low Earth parking orbit. Using the 
generating trajectory as an initial guess, the complexities of the four body problem are added to the targeting 
algorithm incrementally until the final trajectory is converged upon. 

Circular Restricted Three Body Problem 

The reference trajectory for the targeting algorithm is developed in the Sun-Earth circular restricted three body 
problem [16]. In this system, the spacecraft is assumed to be a vehicle of infinitesimally small mass. The motion of 
the spacecraft is governed by two bodies, the Sun and the Earth, both assumed to be point masses. The primary body 
in the system is the Sun, with normalized mass equal to 1-p, and the secondary body is an object with the combined 
mass of the Earth and the Moon located at the Earth-Moon barycenter with normalized mass equal to p, shown in 
Eq. 1. 


Mg + Mm 
Mg + + Mj 


( 1 ) 


The distance unit in the CRTBP is scaled to the average distance from the Sun to the Earth-Moon barycenter. A 
rotating coordinate system is established with the Sun and the Earth-Moon body held fixed on the x-axis, with the 
Earth-Moon body at the origin and the Sun located at x = -1. Assuming the Earth-Moon body orbits the Sun in a 
circular orbit, the coordinate frame rotates at a constant angular velocity, co. The time unit is scaled such that co =1. 
The motion of the spacecraft in the CRTBP defined above is governed by Eq. 2, 
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where the vector r is the position vector in the CRTBP coordinates described above. 

An integral of motion known in the system, the Jacobi integral, is shown in Eq. 3. The mass parameter, p, for the 
system with the Sun as the primary and the Earth and Moon combined into a single secondary body is 
3.040364489 x 10' 6 . 
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Family fl 6 Periodic Orbits 


In the three body problem it is useful to look at periodic orbits as a tool for understanding the dynamics. 
According to the conjecture by Poincare [17], the population of periodic orbits in the CRTBP is dense; ensuring that 
for every orbit there is a periodic orbit an infinitesimally small distance away in phase space. Here, periodic orbits 
around the secondary body (in this system the Earth-Moon combined mass) are used to demonstrate perigee 
increasing effects of the Sun on the spacecraft. Markellos [15] provided a survey of families of periodic orbits in the 
CRTBP. If Markellos’s family /is extended such that the trajectories pass nearer to the secondary body, periodic 
orbits can be generated that have the desired effects on the perigee radii of the orbit. Fig. 1 shows an example of 
such an orbit from family fl6 in Markellos’s catalog. The trajectory is shown in a rotating reference frame centered 
on the Earth-Moon barycenter and with the Sun located on the negative x-axis. It demonstrates the effect of the 
Sun’s gravity on the perigee of an orbit about the secondary mass in the CRTBP. The solar perturbation of the two 
body elliptical orbit increases the perigee radius between the first and second flybys, labeled in Fig. 1 as perigee 1 
and perigee 2, respectively. The magnitude of the third perigee is also greatly increased when compared to the first 
perigee. The periodic orbit has five perigees in total, with perigees 4 and 5 symmetric reflections of perigees 3 and 2, 
respectively. Table 1 compares the perigee radii of the three unique flybys. In Markellos’s scaling of the Sun- Jupiter 
CRTBP, which differs from the scaling in Eq. 1, the orbit has the defining properties listed in Table 2. 




Fig. 1 Periodic orbit from Markellos’s family /16 in rotating coordinates centered on the Earth-Moon 

combined mass 


Table 1 Perigee distances in the orbit /16 

perigee radial distance (km) 

T 7200 

2 230434 

3 132580 

A similar periodic orbit can be found that is nearly a reflection of the orbit in Fig. 1 about the y-axis of the 
rotating coordinate system centered at the Earth-Moon combined mass. The orbit shown in Fig. 2, labeled here /’ 16, 
displays the same perigee raising characteristics as its counterpart in family fl 6. 

An important characteristic of the orbits defined by families fl 6 and/’16 is that members of the family can be 
defined by a single parameter. For example, given a value for the nearest perigee radius that falls in the range of 
existence for the family, a unique member of the family can be defined. The orbits can therefore be scaled so that the 
nearest perigee may coincide with a pre-defined low Earth orbit. Similarly, the Jacobi constant of the orbit can be 
used as a defining parameter. Given a value of the Jacobi energy, a member of each family can be found. 

Table 2 Properties of the fl6 reference periodic orbit in Markellos’s scaling units 

initial x value Jacobi constant 
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Fig. 2 Periodic orbit /’16 in rotating coordinates centered on the Earth-Moon combined mass 


From Miele [18] and Szebehely, [16] an orbit in the CRTBP that makes two perpendicular crossings of the 
x-axis is necessarily periodic due to the symmetries of the system. Both of the orbits detailed cross the x-axis 
perpendicularly at their nearest and farthest points from the secondary body. These periodic orbits are 5-periodic 
orbits of the second kind, crossing the x-axis five times before the periodicity is established with a perpendicular 
intersection with the x-axis. 

Periodic orbits of this type can be easily found using differential correction [19]. The trajectory begins at the 
nearest perigee, located on the x-axis with velocity perpendicular to the x-axis. The trajectory is propagated from 
that point with an initial velocity that yields an appropriate Jacobi constant until the x-axis is crossed five times. The 
initial velocity is then adjusted through an iterative process to yield a trajectory that has zero velocity in the x- 
direction at the 5 th crossing of the x-axis. The converged initial conditions are then propagated forward for twice the 
amount of time, and the periodic orbit is established. 

Reference Trajectory Selection for Ballistic Lunar Capture Targeting 

As detailed in Table 1, this periodic orbit can be used to increase the perigee of a spacecraft’s orbit without a 
maneuver. The reference trajectory for the ballistic lunar capture targeting algorithm will begin with an insertion 
into a family member of the shown periodic orbits at the nearest perigee, labeled perigee 1 in Fig. 1, and end at one 


of the other perigees. The periodic orbit will be selected from family /1 6 or f 16 to coincide at the nearest perigee 
with the radius of the parking orbit of the satellite. Any of the perigees along the periodic orbit that have had their 
radii increased could then be used to target ballistic lunar capture trajectories. The work in this paper, however, 
covers only examples for achieving capture at the first and second perigees. If the first perigee is targeted, the 
transfer trajectory from low Earth orbit (LEO) to lunar orbit will have a period of approximately 100 days, and the 
reference trajectory is shown in bold in Fig. 3. The reference trajectory shown in Fig. 3 is labeled/’/ 6p 1 to denote 
that it is from the generating family /’l 6 and targets a capture at the first increased perigee. In the case of targeting 
the second perigee, the period becomes approximately 180 days, and the reference trajectory, shown in Fig. 4, is 
labeled f’16p2 due to its generating family and its target for capture at the second perigee. It should be noted that 
none of the perigee radii in the periodic orbits are equal to, or even closely approximate, the lunar distance from the 
Earth. These trajectories make suitable reference orbits because they display the dynamics necessary for the perigee 
to increase. Converged solutions that resulted from the reference trajectories of the two different periods are shown 
in Fig. 5 and Fig. 6. The solutions are not geometrically identical to the reference trajectories. The periodic orbits 
simply provide initial guesses that reliably allow convergence to ballistic lunar capture transfers. The selection of the 
appropriate periodic orbit for the reference trajectory from the generating families /1 6 and/’16 will be discussed in a 


following section. 
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Fig. 4 Reference trajectory for a fl6p2 transfer in rotating coordinates centered on the Earth-Moon 

combined mass 
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Fig. 5 A f’16pl transfer in Earth centered Sun-Earth rotating coordinates 
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Fig. 6 A f’16p2 transfer in Earth centered Sun-Earth rotating coordinates 


Targeting Algorithm 


An algorithm to incrementally step from the reference trajectory presented above to a ballistic lunar capture 
trajectory in a real world model is presented here. An incremental approach is used to increase the stability and 
reliability of the algorithm. In each step, elements of the real world system are added, increasing the complexity of 
the dynamics from the CRTBP used in the formulation of the reference trajectory. In this algorithm, the radius of the 
initial parking orbit and the desired initial time of the transfer are given parameters. The converged solution provides 
the orientation of the parking orbit and the AV that will result in a captured lunar orbit. 

Transfer Orientation 

The initial task of the targeting algorithm is to select the proper family of generating obits from the two 
generating families, /1 6 and /I6’. One of the factors of the reliability of the targeting algorithm is the effect of the 
Moon’s gravity on the outbound leg of the ballistic lunar capture trajectory. Depending on the orientation of the 
Moon with respect to the trajectory, the lunar perturbation can vary greatly. To increase the stability of the 
algorithm, a transfer orientation with the farthest distance from the Moon is chosen in order to mitigate the effect of 
the Moon on this leg of the transfer. Trajectories that employ a lunar flyby on the outbound segment can have 
advantages in cost savings, but are not presented here. 

Insertion into either of the two possible transfers can be achieved with a AV in the direction of the velocity 
vector of the spacecraft in a prograde low Earth parking orbit. The selection of the appropriate class of trajectory to 
target on a given launch date is made based on the location of the Moon in the Sun-Earth rotating coordinate frame 
at the time of the transfer. If the Moon is located in the quadrants of the coordinate system farthest from the Sun, 
then the /’ 16 generating family should be used, maximizing the distance between the Moon and the spacecraft. If the 
Moon is located in a quadrant nearest to the Sun, then the generating family fl 6 is appropriate. Fig. 7 and Fig. 8 
illustrate the selection of the reference trajectory. In the example ballistic lunar capture trajectory described below, 
the family /’1 6 is appropriate for the given epoch. 




Fig. 7 Lunar position for generating family fi 6 selection in Earth centered rotating coordinates 



Fig. 8 Lunar position for generating family /’16 selection in Earth centered rotating coordinates 


Estimation of Transfer Time 

The next step in the targeting algorithm is to estimate the transfer time. The initial time is provided among the 


given parameters in the problem statement. The final time is set to coincide with one of the Moon’s crossings of the 


y-z plane in the Earth-Sun rotating coordinates centered on the Earth. Because the Moon’s orbit is out of plane in the 
rotating coordinate system, there are two crossings of the y-z plane each month, one in the positive-y direction and 
one in the negative-y direction. The appropriate crossing is selected to coincide with the reference trajectory. For 
example, in the f’ 16 trajectories the appropriate final time is the time of the Moon’s passage from the third quadrant 
to the fourth quadrant of the rotating coordinate system shown in Fig. 8. In order to allow for the appropriate 
phasing, the final time is selected to be that of the crossing that occurs nearest to 100 days from the given initial 
time. 

Regardless of the initial time, the final time will be chosen based on the orientation of the Moon in rotating 
coordinates. By following this procedure, the initial time is arbitrary, but a spacecraft will arrive at the Moon at one 
of two possible times during each month. 

AV Targeting in the Restricted Three Body Problem 

With an estimation of the transfer time, the first targeting step can be performed. Following the method 
developed by Pu and Edelbaum, [20] an approximate AV 0 for the four body problem is calculated in a three body 
problem by adding the mass of the Moon to the mass of the Earth. The combined mass is located at the Earth-Moon 
barycenter. The radius of the initial parking orbit is then adjusted to produce the same circular velocity around the 
combined mass object as the original parking orbit around the Earth as shown in Eq. 4. 
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This three body problem allows a good approximation of the AV required to insert the spacecraft into the 
transfer. An initial guess of the AV is determined by taking the difference between the spacecraft’s circular velocity 
and its velocity at the closest flyby on the reference trajectory. For this step, the orientation of the parking orbit is 
fixed in the orbital plane of the Earth-Moon barycenter around the Sun, and only the magnitude of the initial burn is 
varied. The AV is assumed to be oriented in the direction of the spacecraft’s velocity vector. 

In the numerical propagation of the trajectory, the Sun and the Earth-Moon combined masses are both treated as 
point masses. Their locations are determined from the JPL DE405 [21] ephemeris. The trajectory of the spacecraft is 
propagated in a non-rotating coordinate frame centered on the Earth-Moon combined mass, and the motion in the 
restricted three body problem (RTBP) is governed by Eq. 5. 
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A differential correction algorithm is used in this targeting step. In the interest of relating this step to following 
steps, the algorithm is defined as follows. Let x p be the vector of parameters that are varied to satisfy the targeting 
constraints. The vector c contains the constraints that must be driven to zero in the numerical routine. The 
differential correction algorithm estimates the derivative of x p with respect to c in order to find the correction of x p 
to make c approach the zero vector. 

For this targeting step, there is only one free parameter, x p = AV. The target is the spacecraft’s crossing of the y- 
z plane of the Earth-Moon barycenter centered, Sun-combined mass rotating coordinate system at tf, c = r x . The 
transversal of the y-z plane will occur in the bottom half of the x-y plane if observed from above, resulting in a 
trajectory similar to the trajectory shown in Fig. 9. 



Fig. 9 Crossing of the y-axis in rotating coordinates centered on the Earth-Moon combined mass 


As can be seen in Fig. 9, the solar influence on the orbit has raised the orbit’s perigee at the second flyby. Table 


3 details the parameters of the first iteration step for the sample trajectory shown in Fig. 9. 


Table 3 First iteration step 

f 0 (Julian date) parking orbit radius (km) A V (km/s) At (days) 

initial 2453611 7200 3.08568 87.42570 

solution 2453611 7200 3.04789 87.42570 


Parking Orbit Orientation Selection in the Restricted Three Body Problem 

With the result obtained above, the approximate AV has been established for the ballistic lunar capture 
trajectory. Continuing with Pu and Edelbaum’s approximation of the four body problem, the approximate 
orientation of the initial parking orbit should be obtained before attempting to propagate the trajectory in the more 
complicated four body dynamical system. 

Let a, P and y be the orientation angles for the parking orbit. These angles define the position of a spacecraft in a 
circular orbit of given radius, and are analogous to the inclination, longitude of the ascending node, and true 
anomaly, respectively. They differ from orbital elements only in that they are defined relative to the x-y plane of the 
solar system as opposed to the equatorial plane of the Earth. In the targeting algorithm, the orientation angles 
represent the orientation of the spacecraft at the time of the transfer insertion burn. The transformation from a set of 
orientation angles to Cartesian coordinates is shown in Eq. 6 . 
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Using the definition of the correction scheme described above, the four free variables in Eq. 7 are iterated from 
their values that were used in the last targeting step to target a point on the y-axis in the Sun-(Earth-Moon) rotating 
coordinate system. The target should be at a distance from the combined Earth-Moon mass equal to the distance 
from the Earth to the collinear Lagrange point of the Earth-Moon three body system exterior to the Moon, 
c = r - r L2 . In this iteration step, four parameters are used to target a three dimensional point. 


x p =(A y a (3 7 ) r 


(7) 



The result of the previous iterations should provide approximate solutions for the orientation of the spacecraft at 
the time of the transfer insertion burn as well as the magnitude of the burn, which is oriented along the direction of 
the velocity vector of the spacecraft. The approximation is satisfactorily accurate to converge to a solution in the 
complex dynamics of the Sun-Earth-Moon four body problem. Table 4 details the results of the second iteration 
step, continuing the targeting of the example trajectory in Fig. 9. 


Table 4 Second iteration step results 



f o (Julian date) 

a (rad) 

P (rad) 

Y (rad) 

AV (km/s) 

At (days) 

initial 

2453611 

0.00000 

2.70878 

0.40907 

3.04789 

87.42570 

solution 

2453611 

0.62247 

1.78275 

0.35308 

3.04889 

87.42570 


Targeting in the Restricted Four Body Problem 

At this stage of the algorithm, the masses of the Moon and the Earth are decoupled, and a more realistic system 
is used. The numerical stability of any targeting algorithm continues to be an issue, and the ballistic lunar capture 
trajectory is targeted in two steps. 

The four free parameters used in the previous targeting step are used again in the four body system. Equations of 
motion governing the four body problem are now used. The three gravitational bodies are treated as point masses 
located at positions determined from the JPL DE405 ephemeris. The coordinate system in which the trajectory is 
propagated is a non-rotating frame centered at the Earth-Moon barycenter. Eq. 8 shows the restricted four body 
problem equations of motion with an indirect term to account for the acceleration of the origin of the coordinate 
system. 
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The target in this iteration is the spacecraft’s radial distance from the Moon. The iteration is considered 
successful if the final radial distance from the Moon is less than the distance from the Moon to L2 of the 
Earth-Moon three body system. In this case the constraint is an inequality constraint, r sc j m < r L ^ m . The parameter 

vector x p is the same as the previous targeting step, Eq. 7. Table 5 shows the details of the four body targeting step. 


Table 5 Restricted four body targeting step 



f o (Julian date) 

a (rad) 

P (rad) 

Y (rad) 

AV (km/s) 

At (days) 

initial 

2453611 

0.62247 

1.78275 

0.35308 

3.04889 

87.42570 

solution 

2453611 

1.31989 

1 .22746 

0.40292 

3.04932 

87.42570 



Energy Minimization in the Restricted Four Body Problem 


The final step of the targeting algorithm is a constrained minimization of the spacecraft’s Keplerian energy with 
respect to the Moon. In this minimization step, let J = KE m be the scalar performance index to be minimized by a 
sequential quadratic programming algorithm [22]. The final time of the transfer is included with the initial time and 
the orientation angles of the parking orbit as free variables, shown in Eq. 9. 


x p =( Av a 7 t f f (9) 

The constraint in this step forces the trajectory to end at a perilune. This constraint, c = /„„ aids in convergence 
to trajectories that do not crash into the Moon and complete an orbit of the Moon without escaping. If desired, an 
additional inequality constraint may be added to ensure the perilune distance is sufficient to avoid collision with the 
lunar surface. 

After a few iterations, the minimization results in a negative Keplerian energy, and thus a captured orbit around 
the Moon. If a particular orbit around the Moon is desired, constraints can be added to the minimization problem 
that will yield the desired orbit, possibly with another small burn. The results of the targeting algorithm are shown in 
Table 6, and the capture trajectory is displayed in Fig. 5. It is noteworthy that the trajectory already possesses 
negative energy with respect to the Moon before this iteration step; however, the minimization of energy is 
necessary to produce a trajectory that remains captured by the Moon for a significant period of time. 


Table 6 Minimization of Keplerian energy 



t 0 (Julian date) 

a (rad) 

P (rad) 

Y (rad) 

EV (km/s) 

A t (days) 

KE (km 2 /s 2 ) 

initial 

2453611 

1.31989 

1.22746 

0.40292 

3.04932 

87.42570 

-0.05379 

solution 

2453611 

1.36268 

1.14152 

0.40272 

3.04971 

99.03939 

-0.10498 


Success Criteria 

A successful targeting algorithm will systematically produce captured trajectories. In a multibody system, 
however, the definition of capture is problematic. In the two body problem, a negative Keplerian energy is sufficient 
for ensuring the spacecraft will not escape the system. Similarly, in the restricted three body problem a Jacobi 
energy value that produces zero velocity surfaces that constrain the motion of the spacecraft to one of the primary 
bodies is an identifier of a captured spacecraft. In the four body problem used in the targeting algorithm, there are no 
constant energy-like quantities that bind the motion of the spacecraft. This problem has been previously dealt with in 



several ways. Belbruno, in defining the weak stability boundary, defines a weakly captured trajectory as a trajectory 
that completes an orbit of a body, returning to the vertical plane from which it began, with negative Keplerian 
energy [1]. In finding WSB transfers using manifolds in the restricted three body problem, Koon et al. define a 
successful transfer as one that transits through the halo orbit into the system of the secondary body [5]. 

The following three metrics for evaluating the capture of a spacecraft into lunar orbit are proposed to evaluate 
the success of the targeting algorithm. At the conclusion of the final step of the algorithm the trajectory is at a 
perilune. Although termination at perilune is not required in the algorithm to have negative Keplerian energy with 
respect to the Moon, a successfully captured trajectory should have negative energy at this point. The first metric 
used to evaluate the algorithm is whether a solution is found that has negative Keplerian energy with respect to the 
Moon at this point. Furthermore, a second, more stringent metric is used that demands a binding of the motion of the 
spacecraft to the Moon. Under the second metric, the algorithm is considered to be successful if there is a second 
perilune with negative Keplerian energy as the trajectory is propagated forward, ensuring an orbit of the Moon. 
Finally, a third metric is used to ensure against collision. Under this metric, the algorithm is deemed to be successful 
if it meets the previous two requirements, and the perilune radii are greater than the Moon’s radius. 

State Transition Matrix Based Derivatives 

The algorithm described above relies on a differential correction algorithm and a sequential quadratic 
programming routine. These methods require the calculation of the derivatives of the target parameters with respect 
to the free variables. Traditionally, finite difference methods have been used in the correction and minimization 
procedures. An alternative approach using analytical gradients computed via the state transition matrix is presented 
here in an attempt to maximize the accuracy of the derivatives and therefore increase the likelihood of convergence. 

Due to the chaotic nature of ballistic lunar capture trajectories, small changes in the initial state can lead to very 
large changes in the final state, making finite difference methods unreliable. The forward difference method has 
been shown to rely on the differencing of several integrations. Errors are introduced through integration error and 
the selection of the numerical step size. If it is too large, truncation errors dominate, while a step size that is too 
small results in round off errors. The appropriate selection of the step size can result in a partial derivative that is 
accurate to one half of the number of digits that the function is known [23]. 



Zimmer and Ocampo [24] have demonstrated the benefits in accuracy from using state transition matrix 
derivatives, as opposed to finite difference methods. First order derivatives relating the states at the initial time to the 
states at the final time can be calculated explicitly by solving for the state transition matrix. The state transition 
matrix relating to the propagation of a spacecraft through a central body force field can be solved for explicitly, and 
its derivatives are exact to first order. The dynamical models used in this algorithm do not yield explicit solutions to 
the state transition matrix, and therefore it must be numerically integrated. Error is introduced in the integration, but 
it is on the order of the error in the integration itself. The result is a more accurate representation of the required 
derivatives. 

The state transition matrix of the three and four body problem is found by numerically integrating the 
relationship in Eq. 11, where F is the time dependent state propagation matrix for the desired dynamical model, 
shown in Eq. 12, where the state vector X is defined in Eq. 13. 


$(t,t 0 ) = F’(t)4 , (t,t 0 ) 


F(t) = 


dX(t) 

dX(t) 


( 11 ) 

( 12 ) 


X = 


(13) 


The state propagation matrix for the three-body and four-body dynamical systems of interest to the algorithm are 
given in Eq. 14, where r is given in Eq. 5 for the three-body dynamical system and in Eq. 8 for the four-body 
dynamical system. 
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For the non-rotating restricted four body problem, the sub-matrix in Eq. 14 that represents the partial derivative of 
the acceleration vector with respect to the radius vector is written as 
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while in the non-rotating restricted three body problem it is 
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Perturbations in the free variables are related to perturbations in the constraint values using the chain rule. The 
targeting algorithm relies on partial derivatives relating both the final state and the final value of the Keplerian 
energy to the five control parameters listed in Eq. 9. For the orientation angle a , the partial derivatives are written in 
the form 


and 


where 


dKE f _ dKE f dX f dX i 
da dXj dX Q da 


dXf _ dX f dX , 
da dX 0 da 


(17) 
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and 
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da 


r (— sin a cos / 3 — cos a sin f3 cos 7 ) 
r (cos a cos j3 — sin a sin (i cos 7 ) 

0 


v (sin a sin (3 — cos a cos j3 cos 7 ) 
v (— cos a sin /3 — sin a cos /? cos 7 ) 
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Similarly, the partial derivatives relating to the orientation angle /3 are written as 


(21) 


dKE f dKE f dX f dX i 
d/3 " dX f dX 0 d/3 


(22) 


and 



(23) 


dX f _ dX f dX i 
d(3 ~ dX 0 d/3 

where the partial derivative of the Keplerian energy with respect to the final state is given in Eq. 22, the partial 
derivative of the final state with respect to the initial state is given as the state transition matrix, and 
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r (— cos a sin [3 — sin a cos / 3 cos 7 ) 
r (— sin a sin (3 + cos a cos / 3 cos 7 ) 
r cos (3 sin 7 
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The partial derivatives relating to the orientation angle 7 are 


(24) 
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and 
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where the partial derivative of the Keplerian energy with respect to the final state is given in Eq. 22, the partial 
derivative of the final state with respect to the initial state is given as the state transition matrix, and 


<9X,. 

cEy 


r sin a sin f3 sin 7 
— r cos a sin (3 sin 7 
r sin f3 cos 7 
v (sin a cos (3 sin 7 ) 
— v (cos a cos f3 sin 7 ) 
v cos [3 cos 7 


The initial maneuver produces the partial derivatives 

dKE f dKE f dX f ax, 
dAV ~ dX f dX 0 dAV 


(27) 
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and 



(29) 


dX f _ dX f dX i 
dAV ~ dX 0 dAV 

where the partial derivative of the Keplerian energy with respect to the final state is given in Eq. 22, the partial 
derivative of the final state with respect to the initial state is given as the state transition matrix, and 

ax, 

dAV 

Finally, the partial derivative of the final state with respect to the final time does not include the state transition 
matrix. 
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and ;i/ is defined in Eq. 8 for the four body problem and Eq. 5 for the three body problem. 
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Results 

Low energy transfers are shown to be found in a systematic way when periodic orbits in the CRTBP are used to 
initialize the algorithm. To validate the robustness of the technique, a computer program was designed to implement 
the algorithm described above with the goal of computing a low energy lunar transfer on a user supplied date. 
Transfers of class fl6p 1 and f’16pl were computed in an automated way. The program was then supplied with a 
randomly generated date and time between January 1, 2010 and January 1, 2012, and a low energy lunar transfer 
was calculated for that date. The program was fed 1000 randomly generated dates and times. The results of the run 
are shown in Table 7. 

Table 7 Algorithm results 

Number of Runs Criteria 1 (negative energy) Criteria 2 (orbit) Criteria 3 (orbit, no collision) 


1000 


100% 


96.50% 


91.10% 




Software Models 


The numerical trajectory propagation was performed assuming the gravitational bodies act as point masses. 
Unless otherwise noted, the locations of the point masses were determined using the DE405 ephemeris. The 
equations of motion as described in the preceding sections are integrated using DLSODA [25]. For both the 
numerical nonlinear equation solving and minimization functions, VF13AD was used. In all simulations, the values 
of the gravitational parameters used are listed in Table 8. 


Table 8 Gravitational parameter values 
Gravitational Parameter Value (km 3 /sec 2 ) 


/G 

/G 

m 


1.32715 x 10 11 
3.986004 x 10 5 
4.9029 x 10 3 


Conclusions 

An efficient, robust method of targeting ballistic lunar capture trajectories is presented. The use of a periodic 
trajectory from family fl 6 in the circular restricted three body problem with the Sun as the primary mass and an 
Earth-Moon combined body as the secondary mass demonstrates the solar effect in raising the perigee of the 
spacecraft’s orbit to approximately lunar radius. This reference trajectory serves as an excellent initial guess in the 
targeting algorithm that produces ballistic lunar capture trajectories on arbitrary launch dates and times. The 
numerical difficulties of targeting chaotic trajectories in the Sun-Earth-Moon four body problems are effectively 
handled by the use of precise gradients obtained through the numerical propagation of the state transition matrix. 
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